home *** CD-ROM | disk | FTP | other *** search
/ EuroCD 3 / EuroCD 3.iso / Programming / vbcc / machines / amiga68k / libsrc / math / math_040 / atan.s < prev    next >
Encoding:
Text File  |  1998-06-24  |  20.6 KB  |  454 lines
  1. *
  2. *   $VER: atan.s 33.1 (22.1.97)
  3. *
  4. *   Calculates the arctangent of the source
  5. *
  6. *   Version history:
  7. *
  8. *   33.1    22.1.97 (c) Motorola
  9. *
  10. *           - snipped from M68060SP sources
  11. *
  12.  
  13.     machine 68040
  14.     fpu     1
  15.  
  16.     XDEF    _atan
  17.     XDEF    @atan
  18.  
  19. *************************************************************************
  20. * atan():  computes the arctangent of a normalized number               *
  21. *                                                                       *
  22. * INPUT *************************************************************** *
  23. *       fp0 = extended precision input                                  *
  24. *                                                                       *
  25. * OUTPUT ************************************************************** *
  26. *       fp0 = arctan(X)                                                 *
  27. *                                                                       *
  28. * ACCURACY and MONOTONICITY ******************************************* *
  29. *       The returned result is within 2 ulps in 64 significant bit,     *
  30. *       i.e. within 0.5001 ulp to 53 bits if the result is subsequently *
  31. *       rounded to double precision. The result is provably monotonic   *
  32. *       in double precision.                                            *
  33. *                                                                       *
  34. * ALGORITHM *********************************************************** *
  35. *       Step 1. If |X| >= 16 or |X| < 1/16, go to Step 5.               *
  36. *                                                                       *
  37. *       Step 2. Let X = sgn * 2**k * 1.xxxxxxxx...x.                    *
  38. *               Note that k = -4, -3,..., or 3.                         *
  39. *               Define F = sgn * 2**k * 1.xxxx1, i.e. the first 5       *
  40. *               significant bits of X with a bit-1 attached at the 6-th *
  41. *               bit position. Define u to be u = (X-F) / (1 + X*F).     *
  42. *                                                                       *
  43. *       Step 3. Approximate arctan(u) by a polynomial poly.             *
  44. *                                                                       *
  45. *       Step 4. Return arctan(F) + poly, arctan(F) is fetched from a    *
  46. *               table of values calculated beforehand. Exit.            *
  47. *                                                                       *
  48. *       Step 5. If |X| >= 16, go to Step 7.                             *
  49. *                                                                       *
  50. *       Step 6. Approximate arctan(X) by an odd polynomial in X. Exit.  *
  51. *                                                                       *
  52. *       Step 7. Define X' = -1/X. Approximate arctan(X') by an odd      *
  53. *               polynomial in X'.                                       *
  54. *               Arctan(X) = sign(X)*Pi/2 + arctan(X'). Exit.            *
  55. *                                                                       *
  56. *************************************************************************
  57.  
  58. ATANA3  dc.l            $BFF6687E,$314987D8
  59. ATANA2  dc.l            $4002AC69,$34A26DB3
  60. ATANA1  dc.l            $BFC2476F,$4E1DA28E
  61.  
  62. ATANB6  dc.l            $3FB34444,$7F876989
  63. ATANB5  dc.l            $BFB744EE,$7FAF45DB
  64. ATANB4  dc.l            $3FBC71C6,$46940220
  65. ATANB3  dc.l            $BFC24924,$921872F9
  66. ATANB2  dc.l            $3FC99999,$99998FA9
  67. ATANB1  dc.l            $BFD55555,$55555555
  68.  
  69. ATANC5  dc.l            $BFB70BF3,$98539E6A
  70. ATANC4  dc.l            $3FBC7187,$962D1D7D
  71. ATANC3  dc.l            $BFC24924,$827107B8
  72. ATANC2  dc.l            $3FC99999,$9996263E
  73. ATANC1  dc.l            $BFD55555,$55555536
  74.  
  75. PPIBY2  dc.l            $3FFF0000,$C90FDAA2,$2168C235,$00000000
  76. NPIBY2  dc.l            $BFFF0000,$C90FDAA2,$2168C235,$00000000
  77.  
  78. PTINY   dc.l            $00010000,$80000000,$00000000,$00000000
  79. NTINY   dc.l            $80010000,$80000000,$00000000,$00000000
  80.  
  81. ATANTBL
  82.         dc.l            $3FFB0000,$83D152C5,$060B7A51,$00000000
  83.         dc.l            $3FFB0000,$8BC85445,$65498B8B,$00000000
  84.         dc.l            $3FFB0000,$93BE4060,$17626B0D,$00000000
  85.         dc.l            $3FFB0000,$9BB3078D,$35AEC202,$00000000
  86.         dc.l            $3FFB0000,$A3A69A52,$5DDCE7DE,$00000000
  87.         dc.l            $3FFB0000,$AB98E943,$62765619,$00000000
  88.         dc.l            $3FFB0000,$B389E502,$F9C59862,$00000000
  89.         dc.l            $3FFB0000,$BB797E43,$6B09E6FB,$00000000
  90.         dc.l            $3FFB0000,$C367A5C7,$39E5F446,$00000000
  91.         dc.l            $3FFB0000,$CB544C61,$CFF7D5C6,$00000000
  92.         dc.l            $3FFB0000,$D33F62F8,$2488533E,$00000000
  93.         dc.l            $3FFB0000,$DB28DA81,$62404C77,$00000000
  94.         dc.l            $3FFB0000,$E310A407,$8AD34F18,$00000000
  95.         dc.l            $3FFB0000,$EAF6B0A8,$188EE1EB,$00000000
  96.         dc.l            $3FFB0000,$F2DAF194,$9DBE79D5,$00000000
  97.         dc.l            $3FFB0000,$FABD5813,$61D47E3E,$00000000
  98.         dc.l            $3FFC0000,$8346AC21,$0959ECC4,$00000000
  99.         dc.l            $3FFC0000,$8B232A08,$304282D8,$00000000
  100.         dc.l            $3FFC0000,$92FB70B8,$D29AE2F9,$00000000
  101.         dc.l            $3FFC0000,$9ACF476F,$5CCD1CB4,$00000000
  102.         dc.l            $3FFC0000,$A29E7630,$4954F23F,$00000000
  103.         dc.l            $3FFC0000,$AA68C5D0,$8AB85230,$00000000
  104.         dc.l            $3FFC0000,$B22DFFFD,$9D539F83,$00000000
  105.         dc.l            $3FFC0000,$B9EDEF45,$3E900EA5,$00000000
  106.         dc.l            $3FFC0000,$C1A85F1C,$C75E3EA5,$00000000
  107.         dc.l            $3FFC0000,$C95D1BE8,$28138DE6,$00000000
  108.         dc.l            $3FFC0000,$D10BF300,$840D2DE4,$00000000
  109.         dc.l            $3FFC0000,$D8B4B2BA,$6BC05E7A,$00000000
  110.         dc.l            $3FFC0000,$E0572A6B,$B42335F6,$00000000
  111.         dc.l            $3FFC0000,$E7F32A70,$EA9CAA8F,$00000000
  112.         dc.l            $3FFC0000,$EF888432,$64ECEFAA,$00000000
  113.         dc.l            $3FFC0000,$F7170A28,$ECC06666,$00000000
  114.         dc.l            $3FFD0000,$812FD288,$332DAD32,$00000000
  115.         dc.l            $3FFD0000,$88A8D1B1,$218E4D64,$00000000
  116.         dc.l            $3FFD0000,$9012AB3F,$23E4AEE8,$00000000
  117.         dc.l            $3FFD0000,$976CC3D4,$11E7F1B9,$00000000
  118.         dc.l            $3FFD0000,$9EB68949,$3889A227,$00000000
  119.         dc.l            $3FFD0000,$A5EF72C3,$4487361B,$00000000
  120.         dc.l            $3FFD0000,$AD1700BA,$F07A7227,$00000000
  121.         dc.l            $3FFD0000,$B42CBCFA,$FD37EFB7,$00000000
  122.         dc.l            $3FFD0000,$BB303A94,$0BA80F89,$00000000
  123.         dc.l            $3FFD0000,$C22115C6,$FCAEBBAF,$00000000
  124.         dc.l            $3FFD0000,$C8FEF3E6,$86331221,$00000000
  125.         dc.l            $3FFD0000,$CFC98330,$B4000C70,$00000000
  126.         dc.l            $3FFD0000,$D6807AA1,$102C5BF9,$00000000
  127.         dc.l            $3FFD0000,$DD2399BC,$31252AA3,$00000000
  128.         dc.l            $3FFD0000,$E3B2A855,$6B8FC517,$00000000
  129.         dc.l            $3FFD0000,$EA2D764F,$64315989,$00000000
  130.         dc.l            $3FFD0000,$F3BF5BF8,$BAD1A21D,$00000000
  131.         dc.l            $3FFE0000,$801CE39E,$0D205C9A,$00000000
  132.         dc.l            $3FFE0000,$8630A2DA,$DA1ED066,$00000000
  133.         dc.l            $3FFE0000,$8C1AD445,$F3E09B8C,$00000000
  134.         dc.l            $3FFE0000,$91DB8F16,$64F350E2,$00000000
  135.         dc.l            $3FFE0000,$97731420,$365E538C,$00000000
  136.         dc.l            $3FFE0000,$9CE1C8E6,$A0B8CDBA,$00000000
  137.         dc.l            $3FFE0000,$A22832DB,$CADAAE09,$00000000
  138.         dc.l            $3FFE0000,$A746F2DD,$B7602294,$00000000
  139.         dc.l            $3FFE0000,$AC3EC0FB,$997DD6A2,$00000000
  140.         dc.l            $3FFE0000,$B110688A,$EBDC6F6A,$00000000
  141.         dc.l            $3FFE0000,$B5BCC490,$59ECC4B0,$00000000
  142.         dc.l            $3FFE0000,$BA44BC7D,$D470782F,$00000000
  143.         dc.l            $3FFE0000,$BEA94144,$FD049AAC,$00000000
  144.         dc.l            $3FFE0000,$C2EB4ABB,$661628B6,$00000000
  145.         dc.l            $3FFE0000,$C70BD54C,$E602EE14,$00000000
  146.         dc.l            $3FFE0000,$CD000549,$ADEC7159,$00000000
  147.         dc.l            $3FFE0000,$D48457D2,$D8EA4EA3,$00000000
  148.         dc.l            $3FFE0000,$DB948DA7,$12DECE3B,$00000000
  149.         dc.l            $3FFE0000,$E23855F9,$69E8096A,$00000000
  150.         dc.l            $3FFE0000,$E8771129,$C4353259,$00000000
  151.         dc.l            $3FFE0000,$EE57C16E,$0D379C0D,$00000000
  152.         dc.l            $3FFE0000,$F3E10211,$A87C3779,$00000000
  153.         dc.l            $3FFE0000,$F919039D,$758B8D41,$00000000
  154.         dc.l            $3FFE0000,$FE058B8F,$64935FB3,$00000000
  155.         dc.l            $3FFF0000,$8155FB49,$7B685D04,$00000000
  156.         dc.l            $3FFF0000,$83889E35,$49D108E1,$00000000
  157.         dc.l            $3FFF0000,$859CFA76,$511D724B,$00000000
  158.         dc.l            $3FFF0000,$87952ECF,$FF8131E7,$00000000
  159.         dc.l            $3FFF0000,$89732FD1,$9557641B,$00000000
  160.         dc.l            $3FFF0000,$8B38CAD1,$01932A35,$00000000
  161.         dc.l            $3FFF0000,$8CE7A8D8,$301EE6B5,$00000000
  162.         dc.l            $3FFF0000,$8F46A39E,$2EAE5281,$00000000
  163.         dc.l            $3FFF0000,$922DA7D7,$91888487,$00000000
  164.         dc.l            $3FFF0000,$94D19FCB,$DEDF5241,$00000000
  165.         dc.l            $3FFF0000,$973AB944,$19D2A08B,$00000000
  166.         dc.l            $3FFF0000,$996FF00E,$08E10B96,$00000000
  167.         dc.l            $3FFF0000,$9B773F95,$12321DA7,$00000000
  168.         dc.l            $3FFF0000,$9D55CC32,$0F935624,$00000000
  169.         dc.l            $3FFF0000,$9F100575,$006CC571,$00000000
  170.         dc.l            $3FFF0000,$A0A9C290,$D97CC06C,$00000000
  171.         dc.l            $3FFF0000,$A22659EB,$EBC0630A,$00000000
  172.         dc.l            $3FFF0000,$A388B4AF,$F6EF0EC9,$00000000
  173.         dc.l            $3FFF0000,$A4D35F10,$61D292C4,$00000000
  174.         dc.l            $3FFF0000,$A60895DC,$FBE3187E,$00000000
  175.         dc.l            $3FFF0000,$A72A51DC,$7367BEAC,$00000000
  176.         dc.l            $3FFF0000,$A83A5153,$0956168F,$00000000
  177.         dc.l            $3FFF0000,$A93A2007,$7539546E,$00000000
  178.         dc.l            $3FFF0000,$AA9E7245,$023B2605,$00000000
  179.         dc.l            $3FFF0000,$AC4C84BA,$6FE4D58F,$00000000
  180.         dc.l            $3FFF0000,$ADCE4A4A,$606B9712,$00000000
  181.         dc.l            $3FFF0000,$AF2A2DCD,$8D263C9C,$00000000
  182.         dc.l            $3FFF0000,$B0656F81,$F22265C7,$00000000
  183.         dc.l            $3FFF0000,$B1846515,$0F71496A,$00000000
  184.         dc.l            $3FFF0000,$B28AAA15,$6F9ADA35,$00000000
  185.         dc.l            $3FFF0000,$B37B44FF,$3766B895,$00000000
  186.         dc.l            $3FFF0000,$B458C3DC,$E9630433,$00000000
  187.         dc.l            $3FFF0000,$B525529D,$562246BD,$00000000
  188.         dc.l            $3FFF0000,$B5E2CCA9,$5F9D88CC,$00000000
  189.         dc.l            $3FFF0000,$B692CADA,$7ACA1ADA,$00000000
  190.         dc.l            $3FFF0000,$B736AEA7,$A6925838,$00000000
  191.         dc.l            $3FFF0000,$B7CFAB28,$7E9F7B36,$00000000
  192.         dc.l            $3FFF0000,$B85ECC66,$CB219835,$00000000
  193.         dc.l            $3FFF0000,$B8E4FD5A,$20A593DA,$00000000
  194.         dc.l            $3FFF0000,$B99F41F6,$4AFF9BB5,$00000000
  195.         dc.l            $3FFF0000,$BA7F1E17,$842BBE7B,$00000000
  196.         dc.l            $3FFF0000,$BB471285,$7637E17D,$00000000
  197.         dc.l            $3FFF0000,$BBFABE8A,$4788DF6F,$00000000
  198.         dc.l            $3FFF0000,$BC9D0FAD,$2B689D79,$00000000
  199.         dc.l            $3FFF0000,$BD306A39,$471ECD86,$00000000
  200.         dc.l            $3FFF0000,$BDB6C731,$856AF18A,$00000000
  201.         dc.l            $3FFF0000,$BE31CAC5,$02E80D70,$00000000
  202.         dc.l            $3FFF0000,$BEA2D55C,$E33194E2,$00000000
  203.         dc.l            $3FFF0000,$BF0B10B7,$C03128F0,$00000000
  204.         dc.l            $3FFF0000,$BF6B7A18,$DACB778D,$00000000
  205.         dc.l            $3FFF0000,$BFC4EA46,$63FA18F6,$00000000
  206.         dc.l            $3FFF0000,$C0181BDE,$8B89A454,$00000000
  207.         dc.l            $3FFF0000,$C065B066,$CFBF6439,$00000000
  208.         dc.l            $3FFF0000,$C0AE345F,$56340AE6,$00000000
  209.         dc.l            $3FFF0000,$C0F22291,$9CB9E6A7,$00000000
  210.  
  211. X       EQU             -12
  212. XDCARE  EQU             X+2
  213. XFRAC   EQU             X+4
  214. XFRACLO EQU             X+8
  215.  
  216. ATANF   EQU             -24
  217. ATANFHI EQU             ATANF+4
  218. ATANFLO EQU             ATANF+8
  219.  
  220. TEMP_SIZE EQU           24
  221.  
  222. ;--ENTRY POINT FOR ATAN(X), HERE X IS FINITE, NON-ZERO, AND NOT NAN'S
  223. _atan
  224.         fmove.d         (4,sp),fp0
  225. @atan
  226.         link            a0,#-TEMP_SIZE
  227.         fmove.x         fp0,(X,a0)
  228.         move.l          (X,a0),d1
  229.         move.l          d1,d0                   ; sign
  230.         move.w          (XFRAC,a0),d1
  231.         and.l           #$7FFFFFFF,d1
  232.  
  233.         cmp.l           #$3FFB8000,d1           ; |X| >= 1/16?
  234.         bge.b           .ATANOK1
  235.         bra.w           .ATANSM
  236.  
  237. .ATANOK1
  238.         cmp.l           #$4002FFFF,d1           ; |X| < 16 ?
  239.         ble.b           .ATANMAIN
  240.         bra.w           .ATANBIG
  241.  
  242. ;--THE MOST LIKELY CASE, |X| IN [1/16, 16). WE USE TABLE TECHNIQUE
  243. ;--THE IDEA IS ATAN(X) = ATAN(F) + ATAN( [X-F] / [1+XF] ).
  244. ;--SO IF F IS CHOSEN TO BE CLOSE TO X AND ATAN(F) IS STORED IN
  245. ;--A TABLE, ALL WE NEED IS TO APPROXIMATE ATAN(U) WHERE
  246. ;--U = (X-F)/(1+XF) IS SMALL (REMEMBER F IS CLOSE TO X). IT IS
  247. ;--TRUE THAT A DIVIDE IS NOW NEEDED, BUT THE APPROXIMATION FOR
  248. ;--ATAN(U) IS A VERY SHORT POLYNOMIAL AND THE INDEXING TO
  249. ;--FETCH F AND SAVING OF REGISTERS CAN BE ALL HIDED UNDER THE
  250. ;--DIVIDE. IN THE END THIS METHOD IS MUCH FASTER THAN A TRADITIONAL
  251. ;--ONE. NOTE ALSO THAT THE TRADITIONAL SCHEME THAT APPROXIMATE
  252. ;--ATAN(X) DIRECTLY WILL NEED TO USE A RATIONAL APPROXIMATION
  253. ;--(DIVISION NEEDED) ANYWAY BECAUSE A POLYNOMIAL APPROXIMATION
  254. ;--WILL INVOLVE A VERY LONG POLYNOMIAL.
  255.  
  256. ;--NOW WE SEE X AS +-2^K * 1.BBBBBBB....B <- 1. + 63 BITS
  257. ;--WE CHOSE F TO BE +-2^K * 1.BBBB1
  258. ;--THAT IS IT MATCHES THE EXPONENT AND FIRST 5 BITS OF X, THE
  259. ;--SIXTH BITS IS SET TO BE 1. SINCE K = -4, -3, ..., 3, THERE
  260. ;--ARE ONLY 8 TIMES 16 = 2^7 = 128 |F|'S. SINCE ATAN(-|F|) IS
  261. ;-- -ATAN(|F|), WE NEED TO STORE ONLY ATAN(|F|).
  262.  
  263. .ATANMAIN
  264.  
  265.         and.l           #$F8000000,(XFRAC,a0)   ; FIRST 5 BITS
  266.         or.l            #$04000000,(XFRAC,a0)   ; SET 6-TH BIT TO 1
  267.         clr.l           (XFRACLO,a0)            ; LOCATION OF X IS NOW F
  268.  
  269.         fmove.x         fp0,fp1                 ; FP1 IS X
  270.         fmul.x          (X,a0),fp1              ; FP1 IS X*F, NOTE THAT X*F > 0
  271.         fsub.x          (X,a0),fp0              ; FP0 IS X-F
  272.         fadd.s          #$3F800000,fp1          ; FP1 IS 1 + X*F
  273.         fdiv.x          fp1,fp0                 ; FP0 IS U = (X-F)/(1+X*F)
  274.  
  275. ;--WHILE THE DIVISION IS TAKING ITS TIME, WE FETCH ATAN(|F|)
  276. ;--CREATE ATAN(F) AND STORE IT IN ATANF, AND
  277. ;--SAVE REGISTERS FP2.
  278.  
  279.         move.l          d2,-(sp)                ; SAVE d2 TEMPORARILY
  280.         move.l          d1,d2                   ; THE EXP AND 16 BITS OF X
  281.         and.l           #$00007800,d1           ; 4 VARYING BITS OF F'S FRACTION
  282.         and.l           #$7FFF0000,d2           ; EXPONENT OF F
  283.         sub.l           #$3FFB0000,d2           ; K+4
  284.         asr.l           #1,d2
  285.         add.l           d2,d1                   ; THE 7 BITS IDENTIFYING F
  286.         asr.l           #7,d1                   ; INDEX INTO TBL OF ATAN(|F|)
  287.         lea             (ATANTBL,pc),a1
  288.         add.l           d1,a1                   ; ADDRESS OF ATAN(|F|)
  289.         move.l          (a1)+,(ATANF,a0)
  290.         move.l          (a1)+,(ATANFHI,a0)
  291.         move.l          (a1)+,(ATANFLO,a0)      ; ATANF IS NOW ATAN(|F|)
  292.         move.l          (X,a0),d1               ; LOAD SIGN AND EXPO. AGAIN
  293.         and.l           #$80000000,d1           ; SIGN(F)
  294.         or.l            d1,(ATANF,a0)           ; ATANF IS NOW SIGN(F)*ATAN(|F|)
  295.         move.l          (sp)+,d2                ; RESTORE d2
  296.  
  297. ;--THAT'S ALL I HAVE TO DO FOR NOW,
  298. ;--BUT ALAS, THE DIVIDE IS STILL CRANKING!
  299.  
  300. ;--U IN FP0, WE ARE NOW READY TO COMPUTE ATAN(U) AS
  301. ;--U + A1*U*V*(A2 + V*(A3 + V)), V = U*U
  302. ;--THE POLYNOMIAL MAY LOOK STRANGE, BUT IS NEVERTHELESS CORRECT.
  303. ;--THE NATURAL FORM IS U + U*V*(A1 + V*(A2 + V*A3))
  304. ;--WHAT WE HAVE HERE IS MERELY  A1 = A3, A2 = A1/A3, A3 = A2/A3.
  305. ;--THE REASON FOR THIS REARRANGEMENT IS TO MAKE THE INDEPENDENT
  306. ;--PARTS A1*U*V AND (A2 + ... STUFF) MORE LOAD-BALANCED
  307.  
  308.         fmove.x         fp2,-(sp)               ; save fp2
  309.  
  310.         fmove.x         fp0,fp1
  311.         fmul.x          fp1,fp1
  312.         fmove.d         (ATANA3,pc),fp2
  313.         fadd.x          fp1,fp2                 ; A3+V
  314.         fmul.x          fp1,fp2                 ; V*(A3+V)
  315.         fmul.x          fp0,fp1                 ; U*V
  316.         fadd.d          (ATANA2,pc),fp2         ; A2+V*(A3+V)
  317.         fmul.d          (ATANA1,pc),fp1         ; A1*U*V
  318.         fmul.x          fp2,fp1                 ; A1*U*V*(A2+V*(A3+V))
  319.         fadd.x          fp1,fp0                 ; ATAN(U), FP1 RELEASED
  320.  
  321.         fmove.x         (sp)+,fp2               ; restore fp2
  322.  
  323.         fadd.x          (ATANF,a0),fp0          ; ATAN(X)
  324.         unlk            a0
  325.         rts
  326.  
  327. .ATANSM
  328. ;--|X| <= 1/16
  329. ;--IF |X| < 2^(-40), RETURN X AS ANSWER. OTHERWISE, APPROXIMATE
  330. ;--ATAN(X) BY X + X*Y*(B1+Y*(B2+Y*(B3+Y*(B4+Y*(B5+Y*B6)))))
  331. ;--WHICH IS X + X*Y*( [B1+Z*(B3+Z*B5)] + [Y*(B2+Z*(B4+Z*B6)] )
  332. ;--WHERE Y = X*X, AND Z = Y*Y.
  333.  
  334.         cmp.l           #$3FD78000,d1
  335.         blt.w           .ATANTINY
  336.  
  337. ;--COMPUTE POLYNOMIAL
  338.         fmovem.x        fp2/fp3,-(sp)           ; save fp2/fp3
  339.  
  340.         fmul.x          fp0,fp0                 ; FPO IS Y = X*X
  341.  
  342.         fmove.x         fp0,fp1
  343.         fmul.x          fp1,fp1                 ; FP1 IS Z = Y*Y
  344.  
  345.         fmove.d         (ATANB6,pc),fp2
  346.         fmove.d         (ATANB5,pc),fp3
  347.  
  348.         fmul.x          fp1,fp2                 ; Z*B6
  349.         fmul.x          fp1,fp3                 ; Z*B5
  350.  
  351.         fadd.d          (ATANB4,pc),fp2         ; B4+Z*B6
  352.         fadd.d          (ATANB3,pc),fp3         ; B3+Z*B5
  353.  
  354.         fmul.x          fp1,fp2                 ; Z*(B4+Z*B6)
  355.         fmul.x          fp3,fp1                 ; Z*(B3+Z*B5)
  356.  
  357.         fadd.d          (ATANB2,pc),fp2         ; B2+Z*(B4+Z*B6)
  358.         fadd.d          (ATANB1,pc),fp1         ; B1+Z*(B3+Z*B5)
  359.  
  360.         fmul.x          fp0,fp2                 ; Y*(B2+Z*(B4+Z*B6))
  361.         fmul.x          (X,a0),fp0              ; X*Y
  362.  
  363.         fadd.x          fp2,fp1                 ; [B1+Z*(B3+Z*B5)]+[Y*(B2+Z*(B4+Z*B6))]
  364.  
  365.         fmul.x          fp1,fp0                 ; X*Y*([B1+Z*(B3+Z*B5)]+[Y*(B2+Z*(B4+Z*B6))])
  366.  
  367.         fmovem.x        (sp)+,fp2/fp3           ; restore fp2/fp3
  368.  
  369.         fadd.x          (X,a0),fp0
  370. .ATANTINY
  371.         unlk            a0
  372.         rts
  373.  
  374. ;--|X| < 2^(-40), ATAN(X) = X
  375.  
  376.  
  377. .ATANBIG
  378. ;--IF |X| > 2^(100), RETURN     SIGN(X)*(PI/2 - TINY). OTHERWISE,
  379. ;--RETURN SIGN(X)*PI/2 + ATAN(-1/X).
  380.         cmp.l           #$40638000,d1
  381.         bgt.w           .ATANHUGE
  382.  
  383. ;--APPROXIMATE ATAN(-1/X) BY
  384. ;--X'+X'*Y*(C1+Y*(C2+Y*(C3+Y*(C4+Y*C5)))), X' = -1/X, Y = X'*X'
  385. ;--THIS CAN BE RE-WRITTEN AS
  386. ;--X'+X'*Y*( [C1+Z*(C3+Z*C5)] + [Y*(C2+Z*C4)] ), Z = Y*Y.
  387.  
  388.         fmovem.x        fp2/fp3,-(sp)           ; save fp2/fp3
  389.  
  390.         fmove.s         #$BF800000,fp1          ; LOAD -1
  391.         fdiv.x          fp0,fp1                 ; FP1 IS -1/X
  392.  
  393. ;--DIVIDE IS STILL CRANKING
  394.  
  395.         fmove.x         fp1,fp0                 ; FP0 IS X'
  396.         fmul.x          fp0,fp0                 ; FP0 IS Y = X'*X'
  397.         fmove.x         fp1,(X,a0)              ; X IS REALLY X'
  398.  
  399.         fmove.x         fp0,fp1
  400.         fmul.x          fp1,fp1                 ; FP1 IS Z = Y*Y
  401.  
  402.         fmove.d         (ATANC5,pc),fp3
  403.         fmove.d         (ATANC4,pc),fp2
  404.  
  405.         fmul.x          fp1,fp3                 ; Z*C5
  406.         fmul.x          fp1,fp2                 ; Z*B4
  407.  
  408.         fadd.d          (ATANC3,pc),fp3         ; C3+Z*C5
  409.         fadd.d          (ATANC2,pc),fp2         ; C2+Z*C4
  410.  
  411.         fmul.x          fp3,fp1                 ; Z*(C3+Z*C5), FP3 RELEASED
  412.         fmul.x          fp0,fp2                 ; Y*(C2+Z*C4)
  413.  
  414.         fadd.d          (ATANC1,pc),fp1         ; C1+Z*(C3+Z*C5)
  415.         fmul.x          (X,a0),fp0              ; X'*Y
  416.  
  417.         fadd.x          fp2,fp1                 ; [Y*(C2+Z*C4)]+[C1+Z*(C3+Z*C5)]
  418.  
  419.         fmul.x          fp1,fp0                 ; X'*Y*([B1+Z*(B3+Z*B5)]
  420.                                                 ;  +[Y*(B2+Z*(B4+Z*B6))])
  421.         fadd.x          (X,a0),fp0
  422.  
  423.         fmovem.x        (sp)+,fp2/fp3           ; restore fp2/fp3
  424.  
  425.         tst.l           d0
  426.         bpl.b           .pos_big
  427.  
  428. .neg_big
  429.         fadd.x          (NPIBY2,pc),fp0
  430.         unlk            a0
  431.         rts
  432.  
  433. .pos_big
  434.         fadd.x          (PPIBY2,pc),fp0
  435.         unlk            a0
  436.         rts
  437.  
  438. .ATANHUGE
  439. ;--RETURN SIGN(X)*(PIBY2 - TINY) = SIGN(X)*PIBY2 - SIGN(X)*TINY
  440.         tst.l           d0
  441.         bpl.b           .pos_huge
  442.  
  443. .neg_huge
  444.         fmove.x         (NPIBY2,pc),fp0
  445.         fadd.x          (PTINY,pc),fp0
  446.         unlk            a0
  447.         rts
  448.  
  449. .pos_huge
  450.         fmove.x         (PPIBY2,pc),fp0
  451.         fadd.x          (NTINY,pc),fp0
  452.         unlk            a0
  453.         rts
  454.